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Based on the spectral statistics obtained in numerical simulations on three dimensional disordered 
systems within the tight-binding approximation, a new superuniversal scaling relation is presented 
that allows us to collapse data for the orthogonal, unitary and symplectic symmetry (f3 = 1,2,4) 
onto a single scaling curve. This relation provides a strong evidence for one-parameter scaling 
existing in these systems which exhibit a second order phase transition. As a result a possible one- 
parameter family of spacing distribution functions, P g (s), is given for each symmetry class /3, where 
g is the dimensionless conductance. 
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The study of critical phenomena is an important sub- 
ject because of the rich variety of systems exhibiting a 
second order phase transition Q. By a second order 
transition we mean a continuous transition between two 
regimes with the the correlation length diverging at the 
transition point. The description of such phenomenon 
leads to the introduction of very important concepts 
such as scaling, renormalization group and universality 
classes. These reflect the fact that the phase transition 
does not depend on the details of the system but only 
on some general symmetries as well as on the dimen- 
sion of the system. A direct consequence is that different 
systems with different Hamiltonians may share the same 
critical exponents, describing the singularity of the phase 
transition, if the symmetry underlying these systems is 
the same and therefore will belong to the same universal- 
ity class. Other features, on the other hand, may be in 
common for different universality classes leaving the pos- 
sibility to derive simple relations between these classes. 
Such a feature is scaling which is exploited in order to 
find the position of the critical point and the value of 
the critical exponent. Even though scaling may be com- 
monplace, the scaling function may be different for the 
different universality classes. 

In this Letter we present a single one-parameter scal- 
ing relation which is common to several different univer- 
sality classes. This relation involves the spectral statistics 
of a three dimensional (3D) disordered system with addi- 
tional degrees of freedom, e.g., strong magnetic field and 
spin-orbit scattering. The choice of this system comes 
from the realization that it exhibits a metal-insulator 
transition (MIT) as a function of the disorder in the ther- 
modynamic limit M. It is generally assumed that the 
critical behavior at the MIT can be classified in terms of 
three different universality classes according to the sym- 
metry of the system: orthogonal [with time reversal sym- 
metry, O(iV)], unitary [without time reversal symmetry, 
e.g. with a magnetic field, U(A^)] and, symplectic [with 



spin orbit-coupling, Sp(A^)]. One then expects different 
critical exponents related to the MIT for the three differ- 
ent universality classes. 

Surprisingly, in spite of the apparent change of uni- 
versality class, the same value of the critical exponent 
has been found, numerically, both in the presence and 
absence of a magnetic field (sB, as well as spin-orbit 
coupling Moreover, Ohtsuki et al. recently showed 

that the anomalous diffusion exponent and also the 
fractal dimension D(2) seem to coincide at the MIT for 
O(N), V(N) or Sp(A^), in agreement with these results. 
It was recently proposed [|6|j8|] that a natural way to un- 
derstand these coincidences would be to invoke the spon- 
taneous breaking of the symmetry right at the MIT. How- 
ever, in a recent paper numerical evidence has been 
presented suggesting a small difference between the scal- 
ing properties of orthogonal and unitary systems. 

The problem is therefore far from being solved and we 
wish to present new evidence concerning how the differ- 
ent universality classes are linked together. This indi- 
cation gives a nontrivial hint about the way in which 
the symmetry parameter enters into the scaling func- 
tion valid for each individual universality class. We also 
present a possible one-parameter family of spacing dis- 
tribution functions, P(s), for each universality class. 

A convenient way to study the MIT is to resort to 
random matrix theory (RMT) and energy level statistics 
(ELS) In RMT the statistics of the energy spec- 

trum are generally described by three different ensembles, 
Gaussian orthogonal (GOE), unitary (GUE), and sym- 
plectic (GSE) depending upon the symmetries mentioned 
above. Recently it has been shown [[^J|,|H],|l3|,[l4| that 
in addition to the two expected statistics, namely either 
GOE, GUE or GSE for the metallic regime and the Pois- 
son ensemble (PE) for the insulating regime, there is a 
third statistics, called the critical ensemble (CE), which 
occurs only exactly at the critical point. 

In order to investigate the MIT we consider the follow- 
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ing tight-binding Hamiltonian Q 

H = ^2e n \n >< n\ + ^V„, m |n >< 



Hi) 



ln(q) +ln(q w ) 



(1) 



with 



V„ 



V 

Vexp(i9 n . m ) 
Vexp(id ntm ) 



orthogonal 

unitary 

symplectic 



(2) 



where the sites n are distributed regularly in 3D space, 
e.g., on a simple cubic lattice. Only nearest neighbor 
interactions are considered. The phase 9 n , m is a scalar 
related to the magnetic field || and 6 n ^ m is a 2 x 2 matrix 
The site energy e„ is described by a stochastic vari- 
able. In the present investigation we use a box distribu- 
tion with variance W 2 /12. The parameter W describes 
the disorder strength and is the critical parameter. 

Based on the above Hamiltonian, the MIT is stud- 
ied by the ELS method, i.e., via the fluctuations of 
the energy spectrum J6t^0|. Starting from Eq. (|l|) 
the energy spectrum was computed by means of the 
Lanczos algorithm for systems of size L x L x L with 
L — 13, 15, 17, 19 and 21 and W ranging from 3 to 
100 averaging over different realizations of the disorder. 
After unfolding the spectra obtained, the fluctuations 
can be appropriately described by means of the spacing 
distribution P(s) |ll[]. This distribution measures the 
level repulsion and is normalized as is its first moment: 
Mi = (s) = 1- 

In order to characterize the shape of P(s), we first cal- 
culate shape descriptive parameters which continuously 
change as we vary external parameters, e.g., the system 
size L or disorder W: 



M 2 



and S str = MS + In Ms, 



(3) 



where fi2 — (s 2 ) is the second moment of P(s), while 
[is = — (sins). These quantities were first introduced to 
describe the spatial-localization properties of general lat- 
tice distributions Jig] and then used for the shape analy- 
sis of P(s) around the MIT |jlq| . It is interesting to note 
that in contrast to previous methods which used only 
part of the information contained in P(s) Jl(],[ll| we con- 
sider here the entire distribution obtained numerically. 
Parameter q is a well-know quantity in probability the- 
ory that describes the peakedness of a distribution func- 
tion. For example, for P(s) = 5(s — 1) we have q = 1. 
The parameter S s t r is called the structural entropy for 
reasons described elsewhere |ll|. These parameters de- 
scribe not only the bulk features of P(s), but also they 
are sensitive to the numerical upper cutoff of the support 
of P(s). 

In order to describe and compare the different univer- 
sality classes within the same method we perform a linear 
rescaling as 



- Hqp) + 

Sstr ~ S\\r 
Sp — Sy/ 



ln(q w ) 
= S 



(4a) 
(4b) 



where index p refers to the PE and w to the Wigner- 
surmise representing the GOE, GUE or GSE respec- 
tively. The choice of such a rescaling defined in Eqs. (H) 
maps the variables S and Q onto the [0, 1] interval, with 
S = Q = 0(1) belonging to the RMT (PE) limit. Fur- 
thermore, in Eq. ([ia| ) it is more natural to use — ln(g) 
instead of q since, similarly to S s tr, it is connected to 
differences of Renyi-entropies JlTj . 
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FIG. 1. Q(L, W) and S(L, W) for the case of symplectic 
symmetry. Continuous curves are polynomial fits. 

As an illustration of the behavior of these parameters, 
in Fig. |l| we report the results for Q(L, W) and S(L, W) 
for the case of spin-orbit coupling (j3 = 4). We can see 
that the data depend on the size of the system except at 
the critical point W c where P(s) is scale invariant. This 
is due to the fact that the MIT is a second order tran- 
sition and that finite-size scaling laws apply close to the 
transition |19[] . These properties were already used with 
success to describe the MIT |l(],|l3|,|o| . In particular it 
was shown that such quantities have a finite size scaling 
behavior and can be written as 

q(L, W) = /(L/foo); S str (L, W) - h(L/^) (5) 

with correlation length £oc(W) ~ \W — W c \~" , and the 
critical exponent v. The functions f(x) and hix) are 
universal in the sense that they do not depend on the 
details of the systems - just on the general symmetries - 
and therefore they directly reflect the universality class 
of the system. 
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From Eq. d|) we can see that, because of the scaling 
behavior of q(L,W) and S s t r {L,W), if we plot S str as 
a function of q we can see the similarities and also the 
differences between the universality classes. The same is 
true for the rescaled parameters Q and S. 

Indeed, Fig. |2| shows clear differences between the or- 
thogonal, unitary and symplectic cases, although all the 
data fall onto special curves irrespective of W and L for 
each case. This figure allows us to determine the scal- 
ing relations for f3 = 1, 2 and 4 without having to derive 
f(L/£,oo), M-k/£oo) and £,oc(W) which are not easy to ob- 
tain numerically due to their singularities at the critical 
point. 
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FIG. 2. S str a function of q for all the symmetry classes. 
(All data are presented in the full range of disorder.) The 
solid curve is obtained from a simple interpolating P(s) due 
to Izrailev Jlq ], The RMT and Poisson distributions appear 
as solid circles. 



sented in Fig. || can be understood with the introduction 
of the dimensionless conductance as a scaling variable. 
We have found that the constant shifts in (|j) for both 
terms — ln(<7) + \n(q\y) and S s t r — Sw, correspond to a 
convolution of different distributions [1231 : 



(6) 



with e x = s. In this case Q — ln(qp) = — ln(gg) — 
ln(gvv) a nd also Sf tr — S^ r + S^ r . 

In Eq. (^J), V g ^(s) is the numerically obtained spacing 
distribution for different symmetry classes parametrized 
by the dimensionless conductance g, which ranges from 
zero to infinity as L and W changes as well, and Wp(t) 
is the RMT limit for g — > oo represented by, e.g., 
the Wigner-surmise. This rescaling provides us with a 
method to study what is beyond the universal level repul- 
sion present in finite systems in the full range of disorder. 
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Using now the rescaling defined above in Eqs. (^), we 
plot S(L, W) as a function of Q(L, W). The results are 
shown in Fig. ||. We see that all the data scale nicely 
onto the same curve indicating the presence of a one- 
parameter superscaling function. The position of the 
MIT moves along the same curve, for f3 — 1,2 and 4, 
as a function of the critical disorder W c which can be 
changed by the magnetic field and spin-orbit scattering 
rate as well as the type of potential scattering. This new 
superscaling relation is very interesting and of impor- 
tance in shedding new light on the MIT in 3D systems. 
New results indicate that the data for (i = 2 in 2D 
scale onto a different curve (see Fig. ||). This point is 
important because it implies that superscaling is not a 
mere consequence of the universality of level repulsion 
but something more subtle. 

Next we will show that the observed relation S(Q) pre- 



FIG. 3. S(L, W) as a function of Q(L, W) for all the sym- 
metry classes. The solid symbols represent the positions of 
the critical points. The continuous (solid, dashed, dotted) 
curves are our analytical estimates, see text for details. For 
comparison the results obtained for the network mode J^J of 
the quantum Hall effect [d = 2, /3 = 2) are also presented. 

The parameters — ln(<7) and S s t r of V g ^{s) give differ- 
ent curves [see Fig. (^)], while after rescaling Q and S 
give the same curve [see Fig. (|)]. This is what is meant 
by the superscaling relation as can be seen in Fig. |^. Scal- 
ing in this context refers to the appearance of g. We will 
show that the parameters Q and S of the function Q g ^{s) 
appearing in Eq. (^) can account for the major part of 
the numerically observed relation. 

In what follows we give an approximate formula for 
'Pg^is) based on analytical calculations and a phe- 
nomenological assumption. First, we point out that 
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Eq. (||) can be solved exactly for the extreme cases of 
a perfect metal (g — > oo) and perfect insulator (g — > 0). 
In the former case the left hand side should equal the 
Wigner surmise and it is easy to show that such a con- 
volution will hold if the Q g .p(s) function as g — ► oo ap- 
proaches a Dirac-delta function, S(s— 1). As for the per- 
fect insulator we have to find Qo,p so that the left hand 
side in each case equals P(s) = cxp(— s). The solution for 
these problems, introducing the notation R(s) = Qo,p(s), 
is [0 



R(y) 



erfc(y) 



/3 = 1 
(3 = 2 



(2y 2 + l)erfc(t/) - ^e"" (3 



(7) 



where y — bs, a — 2/tt, 7r/4, and 97r/64, and b = l/y/w, 
y / 7r/4, and 3y/n/16 for /3 = 1, 2 and 4, respectively. 

These solutions are spacing distributions themselves 
since their zeroth and first moments are normalized to 
unity. The interpolating formula is introduced based on 
the most simple assumption 



Q g (s) 



W{b B 



(8) 



Parameters a g and b g 



are determined from the normal- 
ization conditions (1) = (s) = 1 for each (3. These in- 
terpolating distributions behave in the limit 5 = and 
g — > oo appropriately as defined above. The continuous 



curves in Fig. 
of Q g (s) [Eq 



show that the rescaled Renyi-entropies 
j)] indeed reproduce the results of the 
numerical experiments. Solid, dashed, and dotted lines 
stand for j3 = 1, 2, and 4, respectively. However, we 
see that the analytical curves do not fall onto the same 
curve. This discrepancy may be due to the simplicity of 
the approximation in (|J) and also because of the pres- 
ence of a maximal spacing, i.e. a cut-off in both the 
numerical histogram and consequently in the analytical 
curves. The analytical curves without the upper cut-off 
(not presented here) fall on top of each other within the 
linewidth precision. 

Finally, Fig. || allows us to give an estimate of the 
critical conductance g* . The best fits to the numerical 
histograms give g* = 1.58, 1.46, and 1.34 for j3 = 1, 2, 
and 4, respectively. 

In conclusion, we have presented evidence for a new 
superscaling relation characterizing the MIT in 3D dis- 
ordered systems with different additional degrees of free- 
dom, i.e., in different universality classes. Such a relation 
gives a hint for the derivation of the symmetry depen- 
dence of the scaling function. We have also given an 
approximate analytical formulation of the spacing distri- 
bution where the symmetry parameter (3 and the scaling 
variable g enter in a very clear way. The estimates of the 
critical conductance on the other hand show differences 
for the position of the MIT. This result is complementary 



to the fact that the critical exponent v obtained numer- 
ically in the three cases is the same 

We have to note that in some recent experiments pro- 
viding the same value of the critical exponent [^3|,^4) .as 
well as the absence of influence of the magnetic field |2J] 
and the spin-orbit coupling |2^| at the MIT, show the 
possibility that the superscaling relation presented in this 
Letter could be verified experimentally. 

The method presented in this Letter can be useful in 
the analysis of other phase transitions as well. 
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